↜ Back to index Introduction to Numerical Analysis 2

Solutions for Lecture 3

Exercise 1

In the following I use notation like a_1 to refer to the element a(1) of a Fortran array a.

  1. implicit none
integer i
integer, parameter :: n = 100
real a(n)

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! print the array values to check
print *, a

end
  2. implicit none
integer i
integer, parameter :: n = 100
real a(n)

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! double the values
do i = 1,n
    a(i) = 2 * a(i)
end do

! print the array values to check
print *, a

end

Exercise 2

  1. We want to compute the sum a_1 + a_2 + \cdots + a_n of the elements a_i of the array.

    To compute the sum, we need to perform a sequence of additions to find the partial sums \begin{aligned} s_1 &= 0 + a_1\\ s_2 &= s_1 + a_2\\ s_3 &= s_2 + a_3\\ &\vdots\\ s_n &= s_{n-1} + a_n\\ \end{aligned} Here s_n is the final sum. Since we do not need the intermediate results s_1, s_2, , s_{n-1}, in a program we can reuse the same variable: \begin{aligned} s &\leftarrow 0\\ s &\leftarrow s + a_1\\ s &\leftarrow s + a_2\\ s &\leftarrow s + a_3\\ &\vdots\\ s &\leftarrow s + a_n\\ \end{aligned} Note that in the first step we initialize s to 0. This sequence of operations can be easily expressed by a do loop:

    implicit none
integer i
integer, parameter :: n = 100
real a(n)
real s

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! sum
s = 0.
do i = 1,n
    s = s + a(i)
end do

print *, s

end

    We could also start with s \leftarrow a_1 and perform the do loop only over i = 2, …, n. (Can you think of a reason why the approach in the code above might be better?1)

    Advanced note. Since this is a very common operation, Fortran has a build-in function sum so we could just say s = sum(a).

  2. We already have the sum from the previous program so computing the average is easy by adding the line:

    print *, "average = ", s / n
    implicit none
integer i
integer, parameter :: n = 100
real a(n)
real s

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! sum
s = 0.
do i = 1,n
    s = s + a(i)
end do

print *, "average = ", s / n

end

  3. To have all we need for variance, we only need to add the computation of a_1^2 + a_2^2 + \cdots + a_n^2 (sum of squares) using the same algorithm as in (1.).

    implicit none
integer i
integer, parameter :: n = 100
real a(n)
real s, sq                  ! to store the sum and sum of squares

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! sum and sum of squares
s = 0.
sq = 0.
do i = 1,n
    s = s + a(i)
    sq = sq + a(i) ** 2
end do

! variance
print *, "variance = ", (sq - s**2 / n) / n

end

Exercise 3

  1. As a simple example, let’s consider n =3 first. How do we find the largest value among a_1, a_2 and a_3? In two steps: 1. Find the largest value m_2 among a_1 and a_2; we denote this as \max(a_1, a_2). 2. Find the largest value m_3 among m_2 and a_3, which we denote as \max(m_2, a_3). m_3 is then the largest value among a_1, a_2 and a_3.

    We can continue like this when we have more values a_1, a_2, …, a_n. Therefore to find the largest element of the array with elements a_1, a_2, , a_n, we need to compute “partial” largest values \begin{aligned} m_1 &= a_1\\ m_2 &= \max(m_1, a_2)\\ m_3 &= \max(m_2, a_3)\\ &\vdots\\ m_n &= \max(m_{n-1}, a_n)\\ \end{aligned} Here m_n is the largest element of the full sequence a_1, a_2, …, a_n. As in Exercise 2(1.), we do not need m_1, m_2, …, m_{n-1} and so in a program we can reuse the same variable \begin{aligned} m &\leftarrow a_1\\ m &\leftarrow \max(m, a_2)\\ m &\leftarrow \max(m, a_3)\\ &\vdots\\ m &\leftarrow \max(m, a_n)\\ \end{aligned} Here m \leftarrow \max(m, a_i) means m \leftarrow \begin{cases} a_i, & a_i > m,\\ m, & \text{otherwise}, \qquad \text{(do nothing)} \end{cases} which can be easily implemented using a simple if statement:

    implicit none
integer i
integer, parameter :: n = 100
real a(n)
real m                  ! max

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! assume that the first element is the largest
m = a(1)

do i = 2,n
    ! update m if necessary 
    if (a(i) > m) then
        m = a(i)
    end if
end do

print *, m

end

    Note. The algorithm requires n \geq 1. This is reasonable since the maximum of an empty array is not defined.

  2. To compute the smallest element, we only need to reverse the inequality sign in the condition in 1.

  3. The algorithm is the same as in (1.), we only need to use the absolute values |a_i|, which in Fortran can be computed using abs(a(i)).

    implicit none
integer i
integer, parameter :: n = 100
real a(n)
real m                  ! max

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

m = abs(a(1))
do i = 2,n
    if (abs(a(i)) > m) then
        m = abs(a(i))
    end if
end do

print *, m

end

Exercise 4

  1. To find the index of the first element satisfying a given condition, we need to check the condition for each element a_1, a_2, …, a_n in a sequence, save the first index i for which the condition is satisfied, and stop checking the condition for further elements:

    implicit none
integer, parameter :: n = 100
real a(n)
integer i
integer p                  ! position of element satisfying the condition

do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

p = 0                       ! 0 indicates element not found

do i = 1,n
    if (a(i) > 0) then
        ! we found the element
        p = i
        ! stop the loop
        exit               
    end if
end do

print *, p

end

  2. Since there are no elements in the array that satisfy the condition a_i > 100, we need to report this fact in some way. In the above code, the value 0 is printed. Since the only allowed values of an array index are between 1 and n, we can understand value 0 as “no element satisfying the condition exists in the array”.

Exercise 5

We use the algorithm from Exercise 3 for finding the largest element, but we save the index i whenever we assign to the variable m. The last assignment to m will happen when we see the largest value of an element for the first time.

As an example, consider the array (/ 0, 1, 2, 2, 0 /): the last time we assign a value to m is when we see the first 2. Indeed, the following steps are performed:

  1. Initialization sets m \leftarrow a_1 = 0 (m is updated as so we set p = 1).

  2. Loop i = 2: a_2 = 1 > 0 = m and therefore we assign m \leftarrow a_2 = 1 and save p =2.

  3. Loop i = 3: a_3 = 2 > 1 = m and therefore we assign m \leftarrow a_3 = 2 and save p =3.

  4. Loop i = 4: a_4 = 2 \not> 2 = m and therefore no assignment to m.

  5. Loop i = 5: a_5 = 0 \not> 2 = m and therefore no assignment to m.

In the end, p = 2, the index of the first element with the largest value.

Here’s a Fortran code:

implicit none
integer, parameter :: n = 100
real a(n)
integer i
integer p                   ! position of the first largest element
real m                      ! max

! some initial array
do i = 1,n
    a(i) = min(i**2 - 100, 100 - i)
end do

! first assume that the first element is the largest
m = a(1)
p = 1                       

do i = 2,n
    ! if we find a larger element, we need to update the position of 
    ! the largest element
    if (a(i) > m) then
        m = a(i)
        ! m was assigned to so we save the index i
        p = i
        ! we cannot exit yet since we may still find a larger element
    end if
end do

! index of the first largest element
print *, p

end

Exercise 6

There are many sorting algorithms (ソート).

One of the simplest ones is the selection sort. This algorithm works by observing that if the array a is sorted in ascending order, the value a_i is the smallest value among the values a_i, a_{i+1}, …, a_n for any i= 1, …, n. Therefore to sort the array, we simply iterate for i = 1, …, n - 1 and for each such i find the smallest value a_p among a_i, a_{i+1}, …, a_n and move it to a_i by swapping a_i and a_p. To find p, we use the algorithm from Exercise 5.

Here is an example code:

implicit none
integer i, j, p
integer, parameter :: n = 5
real a(5)
real m

! some initial array of size n
a = (/ 1., 5., -1., 2., -3. /)

! a(i) should be the smallest of the the elements a(i), a(i + 1), ..., a(n)
! for any i = 1, ..., n
do i = 1, n - 1
    ! find the position of the smallest value among a(i), a(i + 1), ..., a(n) 
    p = i
    m = a(i)
    do j = i + 1, n
        if (a(j) < m) then
            m = a(j)
            p = j
        end if
    end do
    ! swap the values a(i) and a(p) = m so that a(i) is the smallest value
    a(p) = a(i)
    a(i) = m
end do

! print the resulting sorted array
print *, a
end

  1. What if the array is empty, i.e., n = 0?↩︎